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Abstract 
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A  calculational  technique  for  use  in  predicting  fallout  far 
downwind  from  nuclear  bursts  is  developed  and  validated.  Possible 
siting  strategies  for  the  next  generation  of  missiles  might  invite  a 
concentrated  attack  by  thousands  of  nuclear  warheads.  Tv  resulting 
fallout  field  could  consist  of  the  superposition  of  thousai  of  single 
burst  patterns.  The  downwind  extent  of  damaging  radiation  /els  would 
extend  beyond  the  distances  to  which  calculations  are  usual  performed 
for  single  bursts.  Numerical  models  currently  available  cannot  be 
extended  to  these  large  downwind  distances  because  of  the  artificial 
pattern  break  up  inherent  in  their  numerical  quadrature  and  because  of 
prohibitive  computing  requirements.  Two  approaches  to  this  problem  are 
taken  here.  First,  a  numerical  smoothing  which  conserves  radioactivity 
is  developed  to  help  prevent  pattern  break  up.  This  is  partially 
successful  in  that  it  extends  the  predictive  range  farther  downwind, 
but  not  far  enough.  The  second  approach  is  to  abandon  the  numerical 
quadrature  —  known  as  disc  tossing  —  and  adopt  a  whole  cloud  smearing 
approach.  The  key  function  needed  for  the  smearing  approach,  the 
fractional  arrival  rate  of  activity  on  the  ground,  is  derived  directly 
from  physical  principles  and  validated  by  comparison  with  an  extensive 

series  of  numerical  (disc  tosser)  predictions. fVThe  prediction  of 

\ 

activity  arrival  rate  appears  to  be  successful;  the  smearing  model 
results  agree  with  the  numerical  method  over  a  wide  range  of  weapon 
yields  and  activity  distributions.  This  function  is  also  employed  in  a 
limited  sensitivity  study  of  particle  size-activity  distribution 


effects  on  ground  arrival  of  activity.  This  sensitivity  study  confirms 
that  the  range  of  carrier  soil  particle  sizes  available  and  the 
distribution  of  bomb  activity  over  those  particles  is  crucial  in 
establishing  the  time  history  of  fallout  arrival.  The  smearing  model 
is  particularly  sensitive  to  fractionation,  defined  here  as  the 
partition  of  activity  between  the  volume  and  the  surface  of  fallout 


particles. 


FAR  FIELD  FALLOUT  PREDICTION  TECHNIQUES 


I.  Overview 


Although  fallout  resulting  from  a  nuclear  attack  against  our 
existing  strategic  systems  would  unquestionably  lead  to  significant 
damage  to  nearby  civilian  populations,  the  impact  would  be  mitigated  by 
the  generally  remote  basing  of  the  strategic  missile  forces. 
Post-attack  fallout  dose  rates  in  the  otherwise  undevastated  regions 
downwind  of  the  burst  point  can  reasonably  be  estimated  by  any  one  of 
several  existing  fallout  models  from  information  about  the  wind  and 
either  yield  or  exposure  surveys.  Possible  future  missile  basing 
schemes,  for  example,  the  hotly  debated  MX  missile  system,  present  less 
tractable  scenarios.  Such  schemes  might  invite  a  concentrated  attack 
by  thousands  of  nuclear  warheads.  The  fallout  field  resulting  from 
such  an  attack  would  consist  of  the  superposition  of  thousands  of 
single  burst  results.  The  downwind  extent  of  lethal  and  incapacitating 
radiation  levels  would  far  exceed  the  distances  for  which  calculations 
are  now  performed  for  single  bursts.  Otherwise  competent  single  burst 
models  now  available  for  use  in  operational  type  studies  cannot  be 
extended  to  these  large  downwind  distances  because  of  the  unrealistic 


fallout  pattern  break-up  inherent  in  their  numerical  quadrature 


schemes,  and  because  of  prohibitive  computer  time  and  memory 


requirements. 

Most  conventional  fallout  models  either  begin  with  or  simulate  the 
formation  of  a  stabilized  debris  cloud  consisting  of  radioactive 
entrained  soil  particles.  The  cloud  is  sectioned  in  the  horizontal 
from  base  to  top,  defining  a  vertical  stack  of  discs  over  ground-zero. 
For  each  disc,  the  constituent  particles  are  grouped  into  a  histogram 
by  size.  One-by-one,  each  group  of  particles  from  each  disc  i3  then 
tracked  through  prevailing  winds  until  it  reaches  the  ground  some 
distance  downwind.  If  the  particle  size  and  cloud  quadratures  are  fine 
enough,  there  will  be  enough  superposition  of  grounded  discs  to  produce 
a  reasonably  accurate  and  smooth  fallout  pattern.  Such  codes  are 
called  disc  tossers. 

DELFIC  is  supreme  among  the  disc  tossers  from  the  standpoint  of 
modeling  competency  (1:7).  In  DELFIC,  the  rise,  expansion,  and 
stabilization  of  the  debris  cloud  and  the  distribution  of  the  entrained 
soil  within  it  are  all  modeled  dynamically.  Cloud  wafers  are  then 
transported  to  the  ground  through  space  and  (potentially)  time-varying 
wind  fields  using  up-to-date  models  for  particle  fall  rates  and 
debris/atmosphere  interactions.  Activity  is  assigned  to  particles  in 
the  grounded  wafers  on  the  basis  of  empirical  studies  of  debris  samples 
from  nuclear  tests.  Contributions  from  specific  nuclides  are  computed 
by  application  of  the  Bateman  equations  to  each  fission  product  mass 
chain  (2:242). 

For  all  its  sophistication,  DELFIC  is  still  plagued  by  one  problem 
inherent  to  all  disc  tossers:  generation  of  a  realistic  fallout 
pattern  requires  quadrature  fine  enough  to  guarantee  that  each  point  of 
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interest  in  the  fallout  field  be  overlaid  by  the  contributions  from 
several  discs.  As  one  follows  the  falling  wafers  farther  and  farther 
downwind,  the  grounding  points  of  different  particle  size  class  wafers 
gradually  separate  until  the  wafers  no  longer  overlap  at  all.  As  a 
result  of  this  phenomenon,  the  MARK-V  version  of  DELFIC  required  about 
10000  wafers  to  adequately  describe  only  the  local  fallout  field 
(3:27).  The  original  attempt  by  DELFIC' s  authors  to  alleviate  this 
problem  consisted  of  replacement  of  the  original  description  of  the 
horizontal  distribution  of  activity  within  each  disc  as  a  constant,  by 
a  gaussian  centered  at  the  cloud  vertical  axis  (1:31).  In  the 
resulting  code,  known  as  the  ASA  version  of  DELFIC,  tops  and  bases  of 
each  monosized  gaussian  cloud  wafer  were  transported  separately  to  the 
ground.  This  eliminated  discontinuities  at  the  wafer  edges  and 
provided  a  greater  horizontal  extent  to  each  wafer,  thereby  postponing 
the  pattern  break-up  problem  to  points  farther  downwind.  As  a  result, 
a  given  prediction  could  be  computed  with  about  half  as  many  wafers  as 
were  required  with  the  MARK-V  version. 

The  approach  taken  by  this  author  to  extend  DELFIC' s  capabilities 
yet  farther  downwind  was  to  use  an  improved  stratagem  for  transporting 
each  wafer,  so  that  its  contents  would  be  smeared  over  a  more  realistic 
area  in  the  fallout  plane.  This  wafer  smearing  or  pattern  smoothing 
technique  is  described  in  Chapter  II.  With  wafer  smearing,  downwind 
voids  were  eliminated,  and  DELFIC' s  quadrature  requirements  were 
reduced  by'  about  a  factor  of  five  over  those  for  the  MARK-V  version  of 
the  code.  Local  fallout  patterns  could  be  calculated  reliably  with 
only  about  2000  cloud  wafers. 


1-3 


Even  with  wafer  smearing,  pattern  distortion  and  break-up  remained 


a  problem  in  attempts  to  calculate  far  field  fallout.  The  logical 
extension  of  wafer  smearing,  whole  cloud  smearing  and  abandonment  of 
the  disc  tossing  concept,  is  described  in  Chapter  III.  In  the 
development  of  this  new  caleulational  technique,  DELFIC  is  discarded 
except  as  a  comparison  standard,  and  variable  winds  and  the  resulting 
curved  pattern  hotlines  are  sacrificed  temporarily  in  the  interest  of 
providing  an  initial  far  field  capability.  In  Chapter  IV,  the  new 
smearing  model  is  validated  by  comparison  of  its  predictions  for  the 
rate  of  activity  deposition  with  DELFIC* s  results  for  several  yields 
and  particle  size-activity  distributions. 

An  advantage  of  the  smearing  model  is  that  it  effectively 


decouples  the  various  input  parameters  that  determine  fallout  from  one 
another.  This  makes  practical  sensitivity  analyses  such  as  that 
presented  in  Chapter  VI. 


II.  Smoothing  PELF I C  Fallout  Patterns 

In  the  ASA  (Atmospheric  Sciences  Associates)  version  of  DELFIC, 
there  is  in  the  cloud  at  stabilization  time  a  wafer  of  debris  particles 
corresponding  to  each  unique  particle  size  class  and  cloud  altitude 
(1:19).  The  vertical  coordinates  of  each  wafer  are  determined 
dynamically  during  the  cloud  rise  simulation  by  the  competing  forces  of 
lift  (produced  by  the  buoyancy  of  the  cloud)  and  of  gravity  on  each 
size  class  of  particles.  The  motion  of  the  smallest  particles  in  each 
size  class  determines  the  stabilization  altitude  of  the  wafer  top;  that 
of  the  largest  particles  determines  the  altitude  of  the  wafer  base. 
After  cloud  stabilization,  however,  the  ASA  version  of  DELFIC  no  longer 
considers  the  particle  size  class  extrema.  Instead,  in  the  advection- 
plus-settling  method  of  transport,  the  geometric  mean  size  for  the  size 
class  is  used  in  the  equations  describing  the  transport  of  the 
size-class/cloud-altitude  wafer  to  the  ground.  The  wafer  top  and  base 
are  transported  separately  (their  initial  altitudes  may  be  considerably 
different),  with  the  mean  particle  size  for  the  size  class  used  in 
computing  the  top  and  base  trajectories.  When  the  wafer  top  and  base 
have  both  impacted,  the  wafer  contents  are  smeared  over  the  area 
between  the  top  and  base  grounding  locations.  The  resulting 
"footprint"  of  an  individual  wafer  conserves  activity  and  reflects  the 
original  vertical  extent  of  the  wafer  and  the  horizontal  distribution 
of  particles  within  the  wafer.  It  fails,  however,  to  retain  any 
information  about  the  extent  in  size  of  the  contributing  particle  size 


class.  Since  this  algorithm  essentially  treats  the  originally 
continuous  particle  size-frequency  distribution  as  a  series  of  Dirac 
delta  functions,  wafers  for  different  size  classes  settle  at  different 
rates.  The  farther  such  wafers  travel  in  the  prevailing  winds  before 
impact  on  the  ground,  the  greater  will  be  the  separation  between  their 
landing  locations.  Voids  eventually  develop  downwind  of  ground-zero 
where  wafers  representing  adjacent  size  classes  fail  to  overlap. 

The  wafer  smearing  algorithm  developed  here  seeks  to  eliminate  the 
downwind  voids  and  thus  to  produce  a  smoother  fallout  pattern  than  is 
possible  with  the  standard  (ASA)  version  of  DELFIC.  The  new  wafer 
smearing  algorithm  is  a  conceptually  simple  extension  of  the  existing 
DELFIC  approach,  although  the  required  coding  modifications  are 
extensive.  The  resulting  code  will  be  referred  to  here  as  the  AFIT 
version  of  DELFIC. 

The  ASA  version  of  DELFIC  traces  the  motion  of  size  class  extrema, 
a  cell  boundary  approach,  only  until  the  end  of  the  cloud  rise 
simulation.  After  that,  it  adopts  a  cell  centered  approach,  tracing 
the  motion  of  wafers  on  the  basis  of  size  class  means.  The  vertical 
extent  of  each  wafer  remains  fixed  until  the  base  impacts  on  the 
ground.  In  contrast,  the  AFIT  model  consistently  employs  a  cell 
boundary  transport  algorithm,  following  the  trajectories  of  size  class 
extrema  from  cloud  stabilization  until  eventual  ground  impact.  Wafer 
tops  are  defined  throughout  transport  by  the  motion  of  the  smallest 
particles  in  the  size  class.  Wafer  bases  are  defined  by  the  motion  of 
the  largest  particles  in  the  size  class.  The  vertical  extent  of  the 
falling  wafer  increases  until  the  base  impacts  on  the  ground,  because 


the  top  and  base  fall  at  different  rates.  Thus,  the  horizontal  extent 
of  the  grounded  wafer  is  defined  by  computing  the  landing  coordinates 
of  the  largest  particles  in  the  size  class  which  originated  at  the  base 
of  a  wafer  and  those  of  the  smallest  particles  in  the  size  class  which 
originated  at  the  top  of  the  wafer.  The  crucial  difference  from  the 
ASA  approach  is  that  wafer  tops  and  bases  are  defined  throughout 
transport  by  the  predicted  motion  of  the  particle  size  extrema  that 
define  the  bounds  of  each  size  class.  Once  a  wafer  has  landed,  the 
ground  distribution  of  its  contents  is  computed  exactly  as  in  the  ASA 
version  of  DELFIC,  by  a  geometric  formalism  involving  grounded  wafer 
base  and  top  parameters  (1:37). 

While  it  might  at  first  appear  that  there  would  be  a  substantial 
penalty  in  this  method  in  terms  of  additional  computer  time  and  memory 
required  to  compute  and  store  additional  wafer  descriptors,  such  is  not 
the  case.  Since  DELFIC  is  an  overlayed  program,  its  memory 
requirements  are  determined  by  the  size  of  the  largest  overlay.  In 
both  the  AFIT  and  ASA  versions  of  DELFIC,  the  largest  overlay  is  the 
output  map  processing  module.  Thus,  the  AFIT  version  requires  no  more 
memory  to  execute  than  does  the  ASA  version.  Moreover,  by  taking 
advantage  of  a  computational  short-cut  which  is  impossible  with  the  ASA 
version,  the  AFIT  version  of  DELFIC  is  able  to  compensate  for  the 
computer  time  required  to  track  its  additional  wafer  descriptors.  In 
transporting  cloud  wafer  tops  and  bases,  the  ASA  version  is  reportedly 
able  to  take  advantage  of  the  vertical  adjacency  of  cloud  wafers  for 
each  particle  size  class,  in  order  to  eliminate  n-1  out  of  what  would 
otherwise  be  2n  transport  calculations  (1:31).  The  code  to  accomplish 


this  simply  does  not,  and  indeed,  cannot  exist  in  the  ASA  version  of 
DELFIC.  Because  the  ASA  version  models  each  particle  cloud  wafer  as  a 
right  circular  cylinder  of  monosized  particles,  the  dual  identities  of 
wafer  tops  and  bases  within  a  single  size  class  are  destroyed  by  wind 
shear  and  stem  wafer  subdivision  before  transport  to  the  ground  is 
begun.  The  AFIT  version  of  DELFIC,  on  the  other  hand,  is  able  to  take 
full  computational  advantage  of  the  dual  identity  that  exists  for  most 
tops  and  bases  of  wafers  belonging  to  ad.lacent  particle  size  classes. 
Except  at  the  top  and  base  of  the  cloud  and  at  the  small  and  large 
particle  extremes  of  the  size  distribution,  and  except  for  wafers 
located  in  the  cloud  stem,  each  wafer  top  and  base  has  a  dual  identity. 
The  description  of  the  top  of  one  wafer  is  identical  to  that  of  the 
base  of  the  wafer  above  it  in  the  cloud  and  belonging  to  the  next 
smaller  size  class.  Since  the  dual  identity  of  all  such  top  and  base 
pairs  is  preserved  by  the  AFIT  version's  cell  boundary  transport 
algorithm,  computing  time  is  saved  by  performing  only  one  calculation 
for  each  unique  wafer  base  and  top.  Preservation  of  this  dual  identity 
that  exists  for  many  of  the  cloud  wafer  tops  and  bases  during  both  the 
drift  that  takes  place  during  cloud  rise  and  throughout  subsequent 
transport  to  the  ground  guarantees  continuity  of  grounded  cloud  wafers. 
Voids  are  eliminated. 

Figure  II— 1  presents  a  comparison  of  dose  rate  contours  from 
fallout  produced  for  a  hypothetical  2.35  KT  surface  burst  by  both  ASA 
and  AFIT  versions  of  DELFIC.  Both  calculations  used  identical  input 
parameters.  The  stabilized  cloud  was  parceled  into  2000  wafers 
consisting  of  100  particle  size  classes  and  20  vertical  cloud  layers,  a 
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Figure  II-1  Fallout  Patterns  for  ASA  and  AFIT  Versions  of  DELFIC 


quadrature  normally  considered  sufficient  for  a  calculation  for  this 
yield  with  the  ASA  version  of  DELFIC.  The  contours  overlay  almost 
perfectly,  except  in  the  immediate  vicinity  of  ground-zero,  where  stem 
fallout  is  significant.  Since  the  radii  of  stem  wafers  depend  on  the 
size  of  the  cloud  cap  at  the  time  of  wafer  separation,  there  are 
generally  great  differences  between  stem  wafer  base  and  top  parameters. 
The  AFIT  version  of  DELFIC  accounts  for  these  differences 
automatically.  In  the  ASA  version,  it  becomes  necessary  to  subdivide 
stem  wafers  in  the  vertical  in  order  to  retain  resolution  of  stem 
diameter  (1:29). 

In  Figure  II-2,  the  ASA  calculation  was  the  same  as  in  the 
previous  contour  map.  The  AFIT  calculation,  however,  used  only  10 
particle  size  classes.  The  patterns  are  still  qualitatively  the  same; 
but  the  lower  level  AFIT  contours  are  shorter  than  their  ASA 
counterparts,  and  the  500  R/hr  contour  is  shifted  east  of  ground-zero. 
This  not  only  demonstrates  that  10  particle  size  classes  is  too  small  a 
number  for  modeling  competency,  but  also  suggests  the  sort  of  behavior 
to  be  expected  at  large  downwind  distances  from  higher  yield  shots. 
Where  there  are  relatively  few  overlapping  grounded  wafers,  the  contour 
patterns  begin  to  reflect  more  the  assumed  horizontal  distribution  of 
activity  in  the  cloud  and  the  formalism  used  to  distribute  that 
activity  on  the  ground,  than  they  do  the  actual  process  of  fallout 
transport  and  deposition. 

Figure  II-3  vividly  illustrates  the  modeling  problem  posed  by 
non-overlap  of  grounded  wafers.  Here,  both  the  ASA  and  AFIT  versions 
used  only  10  particle  size  classes.  Although  contours  for  the  AFIT 
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Figure  II-2  Fallout  Patterns  for  ASA  and  AFIT  Versions  of  DELFIC 


Figure  II-3  Fallout  Patterns  for  ASA  and  AFIT  Versions  of  DELFIC 


version  of  DELFIC  are  distorted  by  this  paucity  of  wafers,  this  example 
shows  that  the  ASA  version  has  lost  all  relationship  to  reality.  The 
structure  in  the  ASA  pattern  is  a  modeling  artifact  produced  by 
discrete  lumps  of  activity  deposited  downwind  of  the  burst  point.  The 
disjoint  segment  of  the  50  R/hr  contour  even  has  the  precise  elliptical 
shape  dictated  by  the  wafer  deposition  algorithm.  These  calculations, 
and  other  examples  not  presented  here,  indicate  that  the  AFIT  version 
of  DELFIC  agrees  with  the  ASA  version  when  both  are  used  with 
sufficiently  fine  quadrature.  Yet,  the  AFIT  version  contours  do  not 
break  up,  as  the  ASA  contours  do,  when  the  quadrature  is  coarse. 

Figures  A-1  through  A-24  (Appendix  A)  present  fallout  contour 
pattern  predictions  of  the  AFIT  version  of  DELFIC  for  a  variety  of 
particle  size-frequency  distributions  for  yields  from  1  to  1000  KT. 
See  Table  II-1  for  details  of  the  distribution  parameter0  Figure  II-4 
shows  DELFIC' s  activity  assignments  for  the  indicated  size-frequency 
distribution  parameters.  All  DELFIC  activity  calculations  assumed 
fission  of  Pu  by  14  MeV  neutrons.  Table  II-2  presents  the  wind 
profile  used.  The  fallout  patterns  are  identified  by  a  code  label 
associated  with  the  chosen  wind  profile,  yield,  and  size  distribution. 
The  wind  information  indicates  date  and  location;  "03W"  means  October 
3,  1968,  at  Washington,  DC.  Yields  are  indicated  in  kilotons,  with  a 
trailing  "S"  for  surface  burst.  The  DELFIC  default  size  distribution 
is  "DEF";  the  WSEG-10,  SIDAC,  or  RAND  is  "WSEG";  the  others  are 
identified  on  the  basis  of  the  particle  size-frequency  median,  for 
example,  "22P7"  for  22.7  Wm.  Thus,  the  identifier,  "03W1000S22P7" , 
refers  to  a  1000  KT  surface  burst  calculated  using  winds  observed  on 


Figure  II-4  Size-Activity  Distributions  for  DELFIC 


Table  II-2  Wash,  DC  Winds  for  3  Oct  1968 


Altitude 

(m) 

Direction 

(deg) 

Speed 

(m/s) 

132. 

174. 

3. 

540. 

206. 

10. 

1049. 

229. 

10. 

1536. 

235. 

11. 

2045. 

233. 

11. 

2570. 

240. 

9. 

3142. 

258. 

11. 

3720. 

265. 

14. 

4376. 

265. 

12. 

5030. 

255. 

10. 

5790. 

255. 

12. 

6580. 

259. 

12. 

7456. 

275. 

15. 

8411. 

275. 

18. 

9482. 

296. 

15. 

10698. 

291. 

24. 

12113- 

285. 

23. 

12954. 

281. 

21. 

13920. 

272. 

17. 

15059. 

268. 

16. 

16440. 

265. 

14. 

17822. 

266. 

9. 

18655. 

249. 

9. 

19623. 

250. 

9. 

20776. 

234. 

9. 

22198. 

250. 

7. 

24039. 

180. 

3. 

25226. 

258. 

3. 

26684. 

271. 

6. 

October  3,  1968,  at  Washington,  DC,  employing  the  size  distribution 
having  a  median  of  22.7  pm  ("High  yield  standard").  Four  of  the  size 
distributions  were  chosen  as  representative  of  a  variety  of  plausible 
distributions  in  a  previous  study  by  Norment  (21:25).  Although 
Norment* s  "High  yield  standard"  was  intended  to  represent  the  RAND  (or 
WSEG-10)  size-activity  distribution  (lognormal  with  a  median  of  89. H  ym 
and  geometric  standard  deviation  of  2.0)  through  the  assumption  that 
the  RAND  data  represented  the  2.5  moment  of  an  underlying  lognormal 
size-frequency  distribution,  the  distribution  identified  as  "WSEG"  in 
Table  II-1  was  developed  here  to  yield  a  DELFIC  activity  assignment 
that  more  closely  mimics  the  RAND  result.  Norment *s  "Low  yield 
standard"  was  chosen  to  represent  the  2.5  moment  of  a  lognormal  fit  to 
the  extensive  particle  size-activity  data  available  from  the  Small  Boy 
test  shot.  The  other  two  distributions  were  similarly  defined  relative 
to  lognormal  size-activity  distributions  in  use  at  the  time  by  the 
Naval  Radiological  Defense  Laboratory.  Relative  to  the  DELFIC  default 
and  "standard"  distributions,  they  represent  rather  extreme 
distributions,  with  large  size-activity  medians,  and  an  enormous  range 
of  particle  sizes. 

While  many  of  the  predictions  appear  quite  reasonable,  it  is 
apparent  that  modeling  artifacts  are  being  introduced  for  some  size 
distributions,  especially  at  low  yields.  Similar  artifacts  would  be 
expected  for  lower  level  contours  (not  shown),  larger  downwind  ranges, 
and  higher  yields.  Although  the  AFIT  version  of  DELFIC  is  able  to 
better  tolerate  a  coarse  quadrature  than  is  the  ASA  version,  it  is 
apparent  that  neither  code  can  make  competent  predictions  at 


arbitrarily  large  downwind  distances.  The  AFIT  version  guarantees  that 
there  will  always  be  some  overlap  between  wafers  which  were  adjacent  in 
the  stabilized  cloud  and  stem;  but  whenever  there  are  only  a  few 
instead  of  many  overlapped  wafers,  the  resulting  contour  pattern  will 
reflect  more  the  assumptions  about  the  horizontal  distribution  of 
activity  in  individual  wafers  than  it  will  a  realistically  expected 
ground  dose  rate.  The  AFIT  version  can  eliminate  the  voids  at  large 
downwind  distances,  but  the  artificial  lumps  and  resultant  structure  in 
contour  patterns  will  remain.  Although  the  range  of  competence  of  any 
disc  tosser  is  somewhat  dependent  on  yield,  wind  structure,  atmospheric 
properties,  and  assumed  particle  size  distribution,  that  range  is 
finite.  The  competence  of  any  disc  tosser  fallout  code  depends  to  an 
extent  on  the  overlap  of  many  cloud  wafers  at  any  point  of  interest  in 
the  fallout  field.  At  distances  so  large  that  it  is  impractical  to 
provide  a  sufficient  number  of  wafers,  some  other  algorithm  must  be 
used  to  predict  the  deposition  of  fallout. 


DELFIC  and  all  other  disc  tosser  fallout  codes  are  limited  in 
ability  to  competently  model  far  field  fallout  deposition  by  the 
requirement  that  many  or,  at  least  several,  cloud  wafers  overlap  r4-  the 
point  or  points  of  interest.  The  ever  finer  particle  size  fil.  ing 
produced  by  minute  differences  in  particle  settling  rates  over  long 
downwind  transit  distances  forces  ever  finer  particle  size  and  cloud 
level  quadrature  on  the  disc  tosser,  if  multiple  wafer  overlap  is  to  be 
insured.  The  stronger  the  wind,  the  more  severe  is  the  overlap 
problem.  Dose  rate  calculations  at  more  than  a  few  representative 
points  far  downwind,  let  alone  any  attempt  to  produce  contour  maps, 
become  expensive  in  terms  of  computer  time,  and  in  the  case  of 
operational  type  studies  requiring  many  thousands  of  calculations, 
impractical.  Smearing  of  discrete,  contiguous  cloud  wafers,  as  is  done 
in  the  AFIT  version  of  DELFIC,  is  at  best  a  partial  solution  to  the 
problem. 

The  prohibitively  fine  particle  size  and  cloud  layer  mesh  called 
for  by  the  far  field  fallout  problem  suggests,  as  a  limit,  that  one 
really  needs  to  smear  the  entire  debris  cloud,  with  its  total  particle 
burden,  in  some  continuous  manner  over  the  downwind  field  —  rather 
than  smearing  individual  wafers.  Although  the  idea  for  a  "smearing" 
fallout  prediction  model  is  not  new,  the  only  such  model  still  in 


widespread  use  today  is  WSEG-10  (4).  First  developed  more  than  20 
years  ago,  it  survives  in  slightly  modified  form  as  the  fallout 
generator  used  by  the  National  Military  Command  System  Support  Center 
in  the  operational  analysis  code,  SIDAC  (Single  Integrated  Damage 
Analysis  Capability)  (5),  and  as  the  fallout  module  for  the  Texas 
Instruments  TI-59  programmable  calculator  nuclear  weapons  effects 
programs  published  by  the  Defense  Nuclear  Agency  (6). 

WSEG-10  probably  owes  its  longevity  more  to  simplicity  and 
calculational  speed  than  it  does  to  modeling  accuracy.  As  a 
non-numeric  model  based  on  empirical  fits  to  old  nuclear  test  data, 
WSEG-10  is  very  fast  and  requires  few  inputs.  Unfortunately,  even  in 
the  absence  of  highly  sheared  or  complex  wind  fields,  WSEG-10  does  not 
agree  all  that  well  with  DELFIC  (7).  DELFIC,  on  the  other  hand,  has 
been  shown  to  be  capable  of  producing  good  agreement  with  most  of  the 
nuclear  test  shots  for  which  sufficient  data  exists  to  provide  the 
required  model  inputs.  Over  the  years,  WSEG-10  has  received  a 
well-deserved  share  of  criticism  by  Russell  (8:209),  Polan  (9:31), 
Norment  (3:103),  and  others,  as  summarized  by  Bridgman  and  Bigelow 
(10).  What  follows  is  not  a  repair  of  the  WSEG-10  model;  rather,  it  is 
a  derivation  of  a  new  fallout  model  based  on  first  principles:  a  model 
capable  of  answering  the  major  criticisms  of  WSEG-10,  and  which  agrees 
with  DELFIC  to  the  extent  permitted  by  modeling  assumptions. 

Smearing  Model  Details 

The  ground  distribution  of  fallout  activity  per  unit  area  at  time 


t,  A(x,y,t),  can  be  expressed  as  that  part  of  the  cloud  activity 
distribution  which  is  deposited  at  (x,y)  from  time  t=0  to  time  t.  If 
A^(t)  is  the  total  activity  in  the  cloud  at  time  t,  then 

A(x,y,t)  s  A  (t)  Jfcf(x,yf (1) 
z  o 

where  f(x,y,f*)  is  the  normalized  horizontal  spatial  distribution  of 
activity  in  the  cloud,  and  g(t"*)  is  the  fractional  arrival  rate  of 
activity  on  the  ground.  Equation  (1)  has  been  called  the  smearing 
equation.  It  is  the  heart  of  any  fallout  smearing  calculation.  If  the 
fractional  arrival  rate,  g(t'),  of  all  cloud  activity  over  all  ground 
area  can  be  found,  and  if  the  horizontal  spatial  distribution  of 
arriving  cloud  activity  per  unit  area,  f(x,y,Ot  is  known,  then  the 
activity,  A(x,y,t),  "smeared"  on  the  ground  at  location  (x,y),  through 
time,  t,  is  the  integral  of  the  product  of  these  two  functions.  By  way 
of  contrast,  a  conventional  disc  tosser  calculation  never  explicitly 
uses  g(t),  but  rather  performs  a  double  numerical  integration  over 
those  particle  sizes  and  initial  particle  locations  which  can 
contribute  to  A(x,y,t).  The  g(t)  function  is  the  key  to  any  smearing 
calculation.  The  WSEG-10  model  proposes  it  is  an  exponential  function 
as  a  result  of  an  empirical  fit  of  some  early  RAND  calculations.  Here 
it  will  be  developed  directly  from  physical  principles. 

Equation  (1),  the  smearing  equation,  usually  is  not  solved  for 
A(x,y,t);  rather,  following  a  series  of  changes  in  dependent  variable, 
one  solves  for  an  artificial  quantity,  the  "unit  time  reference  dose 
rate",  which  is  the  exposure  rate  in  Roentgens  per  hour  or  dose  rate  in 
rads  per  hour  experienced  at  a  detector  3  ft  above  an  infinite  plane 


field  of  uniformly  distributed  debris,  one  hour  after  detonation 
(11:391):* 

Di(x,y)  =  kYF.  fa>f(x,y,t')g(t')dt'  (2) 

i  o 

where  k  is  known  as  the  "source  normalization  constant",  with  units  of 

(exposure/time)  per  (fission  yield/area),  Y  is  the  total  device  yield, 

and  F^.  is  the  fraction  of  that  yield  which  is  due  to  fission,  as 

opposed  to  fusion.  For  Y  expressed  in  kilotons  and  distances  expressed 

in  statute  miles,  consistent  units  for  k  are  (R/hr)/(KT/mi2) . 

Experimental  determination  of  k  on  the  basis  of  observed  fallout 

patterns  is  difficult.  The  choice  of  a  best  value  depends  not  only  on 

physical  considerations,  but  also  on  modeling  assumptions.  For  the 

smearing  model  to  be  derived  below,  however,  the  choice  is  clear. 

Since  the  aim  is  to  reproduce  DELFIC  results  as  closely  as  possible,  k 

is  obtained  directly  from  DELFIC.  Although  k  is  a  weak  function  of 
235  239 

fission  fuel  (  U,  Pu,  etc.)  and  of  the  neutron  energy  spectrum,  a 

mid-range  value  of  2350  (R/hr)/(KT/mi2 )  ,  corresponding  to  fission  of 
239 

Pu  by  14  MeV  neutrons,  may  be  inferred  by  summing  the  activity 
contributions  from  all  DELFIC  particle  size  classes. 

Consider  now  the  horizontal  distribution  of  activity  in  the  cloud, 
f(x,y,t),  the  first  function  in  the  integrand  of  equations  (1)  and  (2). 


The  upper  limit  is  taken  as  infinity,  not  one  hour,  by  convention. 
Thus,  the  reference  dose  rate  should  be  interpreted  as  the  dose  rate  at 
one  hour,  which  would  be  measured  if  all  the  activity  which  would  ever 
be  deposited  at  (x,y)  were  deposited  there  by  one  hour. 


Cloud  geometry  has  seen  various  forms  in  fallout  modeling.  DELFIC 
arrives  at  a  distribution  by  dynamic  modeling  of  the  rise  and 
stabilization  of  the  buoyant,  entraining  fireball.  Most  models  assume 
some  geometric  shape  for  the  cloud  and  determine  dimensions  on  the 
basis  of  empirical  relations  derived  from  test  shot  data.  One  or  two 
right-circular  cylinders  were  early  popular  choices.  WSEG-10  uses 
gaussian  functions  for  both  horizontal  and  vertical  distributions.  The 
WSEG-10  choices  are  reasonable  in  that  they  are  in  agreement  with  what 
little  data  is  available  and  can  be  rationalized  on  the  basis  of  the 
inherently  stochastic  nature  of  atmospheric  transport.  For  this  model, 
the  WSEG-10  functions  are  adopted  to  describe  cloud  geometry  (4:49; 
12).  The  horizontal  distribution  is  given  by 

f (x,y , t)  =  (2n)"*(0 0  I”1  exp{-|( ( (x-Vt) /0  )2+(y/o  ) 2)}  (3) 

x  y  x  y 

where  V  is  an  assumed  constant  effective  wind.  The  initial  vertical 
distribution  is  similarly  taken  as 

h(z)  =  (2tt)”^  0'1  exp{-£( (z-Hc)/0h)2}  (4) 

The  downwind  standard  deviation,  0  ,  is  taken  to  be  a  yield  dependent 
constant  equal  to  the  effective  radius  of  the  cloud,  0o 

0o  =  exp{ .7  +  UnY)/3  -  3.25/(4+UnY+5.4)2)}  (5) 

for  Y  in  megatons  and  0 o  in  statute  miles.  The  crosswind  variance  is 
the  sum  of  the  downwind  variance  and  two  terms  which  account  for 
diffusive  dispersion  and  shear  dispersion: 


o2  =  a|(1+8t/T)  +  (S  a.  t)2 
y  y  n 


(6) 


S  is  the  crosswind  component  of  the  wind  shear.  H  and  a.  (kft),  the 
y  c  n 

mode  and  standard  deviation  of  the  initial  vertical  distribution 
(equation  (4)),  and  T  (hr),  are  yield  dependent  parameters  given  by 

H  =  44.  +  6.UnY  -  .205(£nY  +  2.42)|fcnY  +  2.42 1  (7) 

c  1  1 

a.  =  .  18H  (8) 

n  c 

and 

T  =  (12H  /60  -  2.5(H  /60)2)( 1  -  .5exp{-(H  /25)2})  (9) 

c  c  o 

Although  he  argues  with  the  values  of  H  ,  a0,  and  a.  calculated  by 

c  n 

these  equations,  Norment  (3:80),  in  his  comparison  of  fallout  models, 
concluded  that  the  WSEG-10  treatment  of  diffusive  and  shear  dispersion 
was  adequate.  At  any  rate,  as  the  simple  model  being  derived  here  will 
not  be  capable  of  dealing  with  complex  wind  patterns  as  is  DELFIC;  the 
current  effort  will  focus  on  model  validation  by  comparison  of  fallout 
ground  arrival  rates  with  DELFIC  results.  In  this  phase  of  prediction, 
this  model  seeks  to  be  about  as  good  as  DELFIC. 

As  stated,  g(t),  the  second  function  in  the  integrand  of  equation 
(1)  or  (2),  is  the  keystone  of  any  smearing  model.  It  is  also  the 
source  of  some  of  the  criticisms  of  the  WSEG-10  model  mentioned 
previously.  The  WSEG-10  choice  for  g(t)  is  a  negative  exponential  in 
time  whose  parameters  were  based  on  some  early  RAND  data  (4:4) 


where  T  and  no  are  yield-dependent  empirical  constants  and  K  is 
determined  by  the  normalization  requirement  that  the  integral  of  g(t), 
taken  over  all  time,  be  unity.  With  this  form  for  g(t),  the  gaussian 
form  for  the  horizontal  distribution  of  activity  in  the  cloud,  and  the 
assumption  of  a  constant  effective  wind,  one  conveniently  obtains  an 
integrable  expression  for  the  fallout  "footprint".  Unfortunately,  the 
formulation  also  predicts  sizable  activity  down  within  moments  of  zero 
time.  Norment  has  suggested  that  the  WSEG-10  formulation  for  g(t)  be 
thoroughly  re-evaluated,  with  particular  regard  to  the  parameters  n0 
and  T  (3:78).  Here  g(t)  is  completely  re-derived  without  recourse  to  a 
priori  assumptions  of  functional  form,  activity  distributions,  or  fall 
dynamics. 

Derivation  of  g(t)  and  the  AFIT  Model 

It  is  conventional  in  fallout  calculations  to  model  the  debris 

particles  as  rigid  spheres  falling  at  terminal  velocity  in  a  viscous 

medium,  air.  Although  the  acceleration  of  gravity  is  nearly  constant 

over  altitudes  of  interest,  the  viscosity  of  the  air  varies  with 

altitude  a3  the  atmospheric  density  changes.  As  a  result,  the  particle 

fall  velocity  is  altitude  dependent.  Since  all  the  particles  may 

reasonably  be  assumed  to  have  the  same  density,  the  only  other  variable 

on  which  fall  rate  depends  is  particle  size.  Clearly,  the  altitude  and 

particle  size  dependencies  of  the  fall  velocity  are  separable. 

Therefore,  the  fall  velocity,  V  ,  can  be  replaced  by  the  derivative, 

z 

-dz/dt;  and  the  result  can  be  integrated  over  the  trajectory  of  a 


single  particle  to  obtain  an  expression  for  the  time  required  for  that 
particle  to  impact  with  the  ground: 


V  =  -dz/dt  =  R(r)Z(z) 
z 


(11) 


dz 

Z(z) 


(12) 


Thus,  the  time  of  arrival,  t  ,  is  uniquely  determined  by  the  size 

Sl 

dependence  of  the  fall  dynamics  and  by  the  initial  particle  altitude. 
Said  another  way,  for  each  initial  altitude,  z0,  the  fraction  of 

activity  which  arrives  on  the  ground  within  the  interval  dt  about  t 
must  equal  the  fraction  of  activity  at  that  altitude  which  is 

associated  with  particles  whose  radii  are  in  a  range  dr  about  r,  where 
r  and  t  are  related  by  equation  (12).  That  is,  the  distribution  of 
activity  with  particle  size  is  transformed  by  the  fall  dynamics 
prescription  into  a  distribution  of  activity  with  time  of  arrival: 

|g(t)dt|  =  | A(r)dr|  (13) 

or,  rearranging, 

g(t)  =  A(r) | dr/dt|  (14) 

Thus,  if  A(r),  the  distribution  of  activity  with  particle  size  is  known 
or  assumed,  g(t)  can  be  found,  provided  only  that  | dr/dt |  can  be 

calculated.  Equation  (12)  provides  the  means  to  do  this.  As  it 

stands,  this  equation  is  in  the  form  used  by  conventional  disc  tosser 
codes  to  calculate  the  grounding  times  of  monosized  particle  wafers. 


For  use  in  a  smearing  code,  the  result  is  inverted  to  produce  an 


I* 


expression  for  the  size  of  particles  landing  as  a  function  of  their 
arrival  time.  Symbolically 

r  =  F(t;zo)  (15) 

Differentiation  of  equation  (15)  yields  dr/dt.  For  a  vertically 
distributed  debris  cloud,  equation  (14)  becomes: 

fZT 

g(t)  =  -J  f(z)  A(r)  (dr/dt)  dz  (16) 

ZB 

With  this  result,  the  activity  arrival  rate  depends  explicitly  on  both 
the  activity  distribution  and  the  particle  fall  dynamics  expression, 
functions  which,  in  the  WSEG-10  model,  are  buried  in  empirical  fits  to 
RAND  fallout  data. 


Evaluation  of  equation  (15)  and  its  derivative  requires  selection 
of  a  mathematical  description  of  particle  fall  dynamics.  The  quantity 
of  primary  interest  in  determining  the  selection  is  the  Reynolds 
number,  R,  which  is  the  ratio  of  the  inertial  force  to  the  viscous 
force  experienced  by  a  particle  moving  through  a  fluid  (air): 


R  =  2V  p  r/n  (17) 

z  a 

where  h  is  the  altitude  dependent  dynamic  viscosity  of  air  (kg/m-sec), 
V  is  the  particle  fall  velocity  (m/sec),  r  is  the  particle  radius  (m), 

Z 

3 

and  P  is  the  air  density  (kg/m  ).  When  the  Reynolds  number  is  small 

a 

(say,  R  <  0.05),  viscous  forces  dominate,  and  the  drag  force  on  the 
particle  is  given  by  Stokes'  law: 
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D  =  6irnrV  (18) 

z 

or,  in  terms  of  the  Reynolds  number,  from  equation  (18) 

D  =  37rn2R/p  (19) 

a 

For  a  particle  falling  at  terminal  velocity  (steady  state),  the  drag 
force  is  balanced  by  the  gravitational  force,  so  that,  for  particles  of 
density  (taken  here  to  be  2600  kg/m3) 

6TrnrVz  =  (4/3)nr3gPf  (20) 

where  g  is  the  gravitational  acceleration  (9.8  m/sec2),  and  buoyancy  is 
neglected,  since  p  »  p  .  An  approximate  cut-off  particle  size  for 

I  3 

which  Stokes'  law  is  valid  is  determined  by  using  equation  (19)  to 
express  equation  (20)  in  terms  of  the  Reynolds  number,  and  solving  for 
the  particle  size: 


r 


=  0 


9h2R  > 
4gPfPaj 


1/3 


(21) 


The  cut-off  size  is  about  10  pm  with  an  admittedly  pessimistic  0.05  for 
the  Reynolds  number,  and  typical  sea  level  atmospheric  properties.  As 
particles  ten  to  a  thousand  times  larger  contribute  significantly  to 
local  fallout,  it  becomes  necessary  to  consider  the  regime  of  larger 
Reynolds  numbers,  where  dynamic  pressure  becomes  dominant  in 
determining  the  drag  force.  For  spheres  the  Stokes  expression  becomes 


D  =  (1/2)p  V2C  .TTr2 
a  z  d 


(22) 
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2 

where  (1/2)p  V  C  is  the  dynamic  pressure  exerted  on  the  projection  of 

cl  2  Q 

the  particle  area  in  the  direction  of  motion  and  is  the  drag 
coefficient.  In  terms  of  R 

D  =  irn2CdR2/(8pa)  (23) 

The  large  Reynolds  number  analog  of  equation  (20)  becomes 

(1/2)PaV2Cdur2  =  (H/3)Trr3gpf  (24) 

or,  using  equation  (23)  in  lieu  of  (22),  and  solving  for  the 

dimensionless  R2C. 

d 

R2Cd  =  32gPaPfr3/(3h2)  (25) 

DELFIC  calculates  V  on  the  basis  of  equations  (17)  and  (25)  by  a 

z 

method  due  to  McDonald  (13)  and  Davies  (14).  R2C.  is  calculated  from 

a 

equation  (25).  The  Reynolds  number  is  then  obtained  from  a  polynomial 
expansion  in  R2Cd  suggested  by  Davies.  Finally,  the  fall  rate  is 
obtained  from  equation  (17).  In  the  interest  of  preserving  maximum 
equivalence  between  the  AFIT  smearing  model  and  DELFIC,  the  same 
approach  is  taken  here.  The  difficulty  with  this  method  is  that  it 
does  not  permit  an  analytic  solution  to  equation  (15).  Thus,  in  order 
to  better  illuminate  the  process,  and  to  indicate  the  small  particle 
limit  which  the  more  exact  solution  must  approach,  it  is  instructive  to 
consider  first  the  solution  for  particles  which  obey  Stokes'  law, 
equation  (20). 

When  Vz  is  replaced  by  -dz/dt  in  equation  (20)  and  the  result  is 
integrated  over  the  trajectory  of  a  single  particle,  the  following 


expression  relating  r  and  t  is  obtained: 


9  Jzn(z')dz' 

2  _  _ o _ 

2gPft 


(26) 


Differentiation  with  respect  to  time  at  constant  altitude  yields 
immediately 


dr/dt  =  -r/(2t)  (27) 

which  can  be  expressed  in  terms  of  t  alone  by  means  of  equation  (26). 
Of  course,  equation  (26)  can  also  be  used  to  express  the  activity 
distribution,  A(r),  as  a  function  of  t  in  lieu  of  r.  Both  of  these  "t 
for  r"  substitutions  will  be  made  later  in  calculating  g(t)  with  use  of 
the  McDonald-Davies  method. 

Outside  of  the  small  particle  Stokes'  law  regime,  equation  (24) 
cannot  simply  be  integrated,  as  was  equation  (20),  to  obtain  a  high 
Reynolds  number  analog  of  equations  (26)  and  (27).  This  is  because  the 
drag  coefficient  is  itself  a  velocity  dependent  variable.  The  method 
of  McDonald  and  Davies,  however,  provides  a  means  to  calculate  fall 
velocities  indirectly  through  the  Reynolds  number,  R.  The  quantity 
R2Cd  is  obtained  from  equation  (25)  and  used  to  obtain  R  by  use  of 
Davies'  polynomials: 

R  =  R2C  , /24  -  2.3363x10"1<(R2Cj)2  +  2.0154x10~6 (R2C  .) 3 
d  q  d  d 

-  6.9105x10"y(R2C.)"  for  R  <  4  and  R2C  <  140 

a  a 

logR  =  -1.29536  +  0.986llog(R2Cd)  -  0.046677(log(R2Cd) )2 

+  0.0011235(log(R2Cd)]3  for  3  <  H  <  104  (28) 
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Note  that  in  the  small  R  limit,  this  formulation  approaches 
R  =  24/C^  ,  which  reduces  the  Newtonian  drag  expression,  equation  (23), 
to  Stokes'  law,  equation  (19).  The  Reynolds  number  thus  obtained  is 
used  in  equation  (17)  to  obtain  the  particle  fall  rate.  McDonald 

originally  suggested  this  algorithm  with  reference  to  a  graphical 

2  2 
representation  of  R  C .  ,  and  observed  that  the  formulation  for  R  C . 

a  a 

given  by  equation  (25)  could  be  readily  exploited  to  calculate  altitude 

dependent  fall  velocities,  since  in  it,  R2C^  depends  only  on  material 

properties  of  the  particle  and  of  the  atmosphere.  DELFIC's  authors 

adopted  the  Davies'  polynomials  in  lieu  of  the  graphs.  The  fall  rates, 

V  ,  thus  calculated  are  corrected  for  high  altitude  drag  slip  by 
z 

multiplying  them  by  the  same  factor  as  is  used  in  DELFIC  (15:6): 

Cf  =  1  +  1.l65x10"7/(Par)  (29) 

for  r  and  P  in  MKS  units.  The  fall  velocities  calculated  in  this  way 
a 

for  each  particle  size  can  be  used  to  numerically  integrate  equation 
(24)  over  the  particle  trajectory,  leading  to  a  numeric  expression  of 
equation  (15),  which  is  analogous  to  the  Stokes'  law  expression, 
equation  (26). 

Since  there  is  no  reason  to  repeat  the  above  procedure  for  each 
fallout  prediction,  data  satisfying  equation  (15),  and  which  span  the 
parameter  space,  are  generated  once.  Th„  .^suiting  r  =  F(t;z<>) 
curves  are  then  fit  to  a  Laurent  series  in  time.  Since  the  Stokes 

results  above  show  that  in  the  small  particle  size  limit  (long  fall 

-1/2  -3/2 

time),  r  is  proportional  to  t  and  dr/dt  is  proportional  to  t  , 


correct  behavior  in  this  limit  is  guaranteed  by  using 

r  =  f  C. (z)ti_^  +  C?(z)t~^  (30) 

i=1  ' 

The  derivative,  dr/dt,  is  found  trivially  by  a  term-by-term 
differentiation  of  the  series.  The  complete  table  of  coefficients, 
spanning  the  0  to  50  km  altitude  range  in  intervals  of  200  m  is 
presented  in  Table  B-1  (Appendix  B).  An  identical  set  of  coefficients 
was  previously  reported  by  Colarco  (16:67). 

The  calculation  of  A(r) dr/dt  in  the  integrand  of  equation  (16) 
proceeds  as  follows  for  each  altitude  required: 

1 .  Interpolate  in  the  table  of  coefficients  for  the 
Laurent  series  expansion  of  r  =  F(t;z0)  to  obtain  coefficients  for 
the  current  altitude,  z0. 

2.  Evaluate  the  series  expansions  for  r  and  dr/dt  at 

time  t. 

3.  Use  the  calculated  value  of  r  to  evaluate  A(r). 

The  limits  of  integration  for  equation  (16)  are  established  as 

±3cr  about  the  cloud  center  height,  H  ,  as  given  by  equation  (7).  The 
n  c 

vertical  distribution  of  activity  is  a  gaussian,  with  standard 
deviation,  =  0.l8Hc  (equation  (8)). 

Figure  III-1,  which  presents  g(t)  for  a  1000  KT  yield,  compares 
the  forms  of  g(t)  predicted  by  WSEG's  negative  exponential,  a  RAND 
model  calculation  reported  in  RM-2460  (17),  the  AFIT  smearing  model 
developed  above  (all  three  using  the  same  RAND  lognormal  size-activity 
distribution).  In  addition,  a  DELFIC  calculation  is  presented,  which 
uses  an  approximation  to  the  RAND  activity  distribution  (introduced  in 
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Figure  III-1  Deposition  Rates  for  RAND's  Activity  Distribution 
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the  previous  chapter).  The  WSEG  and  AFIT  models  used  the  same  cloud 
model;  RAND's  cloud  was  a  right  circular  cylinder.  Although  there  are 
minor  differences  in  magnitudes,  the  RAND,  DELFIC,  and  AFIT  models  all 
predict  the  same  general  behavior  for  the  g(t)  function.  The  two  WSEG 
calculations  are  notably  different,  both  at  early  and  late  times.  Both 
WSEG  calculations  used  a  decay  time  of  7.5  hours  for  the  activity 
deposition  rate,  as  determined  by  the  SIDAC  prescription  (5:14).  The 
no  values  of  1  and  2  bracket  the  range  of  possible  values  for  all 
yields.  WSEG's  choice  of  a  functional  form  for  its  g(t)  fit  leads  to 
the  prediction  of  sizable  activity  down  within  moments  of  zero  time. 
This  early  over-prediction  is  compensated  for  by  a  too  rapid  decay  in 
the  activity  grounding  rate  at  later  times.  Since  the  size-activity 
distributions  and  cloud  geometries  were  the  same  for  both  WSEG  and  AFIT 
models,  the  AFIT  model  apparently  offers  an  improvement  over  WSEG  in 
modeling  fallout  arrival  rate. 

Although  most  fallout  models  use  a  lognormal  function  to  depict 
the  distribution  of  activity  with  particle  size,  as  was  done  for  the 
comparison  above,  there  is  no  need  for  such  a  restriction  with  the  AFIT 
model.  Part  of  DELFIC' s  flexibility  derives  from  its  avoidance  of  such 
an  a  priori  selection.  Rather,  on  default  of  input  to  the  contrary, 
DELFIC  specifies  a  lognormal  frequency  distribution  of  particle  size, 
in  accord  with  Freiling's  extensive  observations  of  fallout  particle 
sizes  (18:6): 

dN(r)  =  dr -  exp{-i(^— TT"*)  }  (3D 

/2?r  ins  *ns 
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where  r0  is  the  median  radius  and  s  is  the  geometric  standard  deviation 
(dimensionless).  This  form  is  especially  convenient  because  of  the 
simple  relationships  which  exist  between  the  various  moments  of  a 
lognormal  distribution  (19:12).  Specifically,  the  n-th  moment  is  also 
lognormal,  with  the  same  geometric  deviation  and  with  a  median  given  by 

rn  =  exp{S-nro  +  n(i-ns)2}  (32) 

Thus,  surface  area  is  lognormally  distributed  with  Jlnr,,  =  &nr0  + 
2 (Ins) 2  ,  and  volume  is  distributed  with  Jtnr^  =  fi-nr0  +  3(#>ns) 2  . 
DELFIC  then  models  the  assignment  of  activity  to  the  particle 
distribution  according  to  a  modification  of  Freiling’s  radial  power 
distribution  law,  which  assumes  that  the  nuclides  of  each  fission 
product  mass  chain  are  distributed  according  to  some  non-integer  moment 
of  the  size-frequency  distribution,  some  of  each  chain  being  mixed 
homogeneously  throughout  the  carrier  soil  particle  (the  "refractory" 
fraction),  the  remainder  being  plated  out  on  the  surface. 
Fractionation  is  the  term  generally  used  to  describe  this  process.  The 
aggregate  size-activity  distribution  is  just  the  sum  over  all  mass 
chains. 

In  this  model,  DELFIC’ s  use  of  a  particle  size-frequency 
distribution,  is  retained,  but  the  individual  mass  chain  detail  is 
omitted.  Even  so,  the  following  simple  assumption  permits  inclusion  of 
a  first  order  accounting  for  fractionation  which  is  based  on  DELFIC 
results.  A  fraction,  fp,  of  the  total  cloud  activity  is  assumed  to  be 
distributed  homogeneously  throughout  the  carrier  soil  particles,  or 
according  to  the  third  moment  of  the  size  frequency  distribution.  The 


remaining  fraction,  (1-fp),  is  assumed  to  be  plated  out  on  the  surface 
of  the  carrier  soil;  and,  hence,  it  is  distributed  according  to  the 
second  moment.  The  total  activity  distribution  is  then  the  sum  of 
volumetric  and  surface  components: 

dA(r)  .  (f  expl-it*"--  - 

/2rir  Jins  r  *ns 

-h  (1  -  fr)  exp{-i(^-^^)2})  (33) 

DELFIC's  assignments  of  activity  to  particle  size  groups  provide  a 

means  to  estimate  f^.  Equation  (3*0  is  integrated  to  produce  an 

expression  for  f  in  terms  of  the  cumulative  3ize-activity 

distribution,  Ac(r),  and  the  cumulative  surface  and  volumetric 

components,  A  (r),  and  A  (r),  respectively: 
vc  sc 

fr  =  (»c<r)  -  AS0<D)  /  (Avo(r)  -  As(,(r))  <34> 

where  corresponding  values  of  A  (r)  and  r  are  obtained  from  DELFIC. 

c 

Since,  as  noted  above,  DELFIC  does  not  use  a  simple  two-component 

distribution,  values  of  fp  calculated  in  this  way  are  weakly  dependent 

on  the  choice  of  A  (r)  and  r.  Selection  of  the  median  radius  of  the 

c 

DELFIC  size-activity  distribution,  the  size  for  which  AQ(r)  =  0.5  , 
yields  a  two-component  distribution  which  closely  follows  DELFIC  over 
the  entire  particle  size  range.  Resulting  f  values  (see  Figure  III-2) 
were  found  to  be  essentially  independent  of  yield  over  a  1  KT  to 
1000  KT  yield  range,  and  only  weakly  dependent  on  choice  of  particle 
size-frequency  distribution.  Figures  III-3  through  III-8  display  the 
two-component  distributions  derived  in  this  manner  for  each  of  the 


Figure  III-2  Refractory  Activity  Fractions 
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size-frequency  distributions  introduced  in  Chapter  II.  Also  shown  are 
the  contributing  surface  and  volumetric  distributions  and  the  DELFIC 
data  on  which  the  calculations  of  f  for  the  two-component  distribution 
were  based. 

Approximations  to  the  Basic  AFIT  Model 

In  the  preceeding  derivation  of  the  AFIT  model,  g(t)  is  expressed 
as  an  integral  over  the  vertical  cloud  distribution.  Unfortunately, 
this  integral  can  only  be  evaluated  numerically.  The  model  may  be 
further  simplified,  however,  by  a  suitable  approximation  of  the 
integration  over  the  vertical  cloud  extent. 

If  the  vertical  cloud  distribution  is  replaced  by  a  Dirac  function 
at  Hc>  the  cloud  degenerates  into  an  infinitesimal  wafer  at  the  center 
height,  and  the  integration  becomes  a  simple  matter  of  evaluating  the 
integrand  at  Hc«  The  activity  arrival  rate,  g(t),  is  then  given  by 
equation  (14),  where  both  A(r)  and  dr/dt  are  evaluated  as  described 
above,  with  z  =  Hc.  Figure  III-9  compares  AFIT  model  g(t) 
calculations  which  used  both  this  "pancake"  and  vertically  distributed 
gaussian  clouds  for  1,  10,  100,  and  1000  KT  yields  and  the  DELFIC 
default  size-frequency  distribution.  Colarco  (16:87)  presented  similar 
results  for  lognormal  size-activity  distributions  and  other  yields 
ranging  from  1  to  20000  KT.  In  all  cases,  there  are  negligible 
differences  between  "pancake"  cloud  and  distributed  cloud  model 
results.  This  occurs  because  the  integrand  in  equation  (16)  is  nearly 
symmetric  between  the  integration  limits  about  Hc>  It  should  be 
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remembered,  however,  that  any  prediction  based  on  a  stabilized  cloud 
cap  model  such  as  AFIT  or  WSEG,  can  be  expected  to  give  unrealistic 
predictions  close  to  ground  zero  for  a  low  yield  explosion,  because  of 
neglect  of  stem  fallout.  It  is  also  important  to  note  that  the 
"pancake”  approximation  is  used  only  with  respect  to  evaluation  of  the 
integral  over  altitude  in  the  g(t)  expression.  Fallout  footprint  or 
dose  calculations  that  might  be  based  on  g(t)  calculated  in  this  manner 
would  still  use  a  distributed  cloud  in  modeling  the  effects  of 
horizontal  wind  shear;  specifically,  the  vertical  standard  deviation  of 
the  cloud  would  still  taken  as  0.1 8Hg*  In  the  following  chapter,  this 
"pancake"  approximation  is  used  in  validating  the  AFIT  model 
predictions  for  g(t)  against  DELFIC's  results. 

Regardless  of  the  method  used  to  obtain  g(t),  whether  by  numerical 
integration  or  by  use  of  the  "pancake"  approximation,  equation  (2)  for 
the  unit  time  reference  dose  rate,  in  which  g(t)  is  used  to  smear  the 
cloud  distribution  over  the  ground,  involves  an  integration  of  g(t).  A 
suitable  analytic  approximation  for  this  integration  is  obtained  by 
expansion  of  g(t)  ir.  a  Taylor's  series  about  the  downwind  time  of 
arrival,  t  ,  of  the  cloud  center  at  location  x.  Retention  of  only  the 
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first  two  terms  yields 
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a  a 

(t) 

Ift 

a 

(35) 

When  this 

approximation 

is 

used  in  equation  (2),  and 

the  change  of 

variable, 

u  =  -(x-Vt)/o 

9 

is  made,  the  integral  in 

that  equation 

becomes 
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x 

where  f(u)  and  f(y,t  )  represent,  respectively,  the  downwind  and 
crosswind  components  of  f(x,y,t).  Now,  in  the  domain  x  =  Vt  >30  , 

cl  X 

cloud  passage  begins  after  zero  time,  and  virtually  no  error  is 
introduced  by  extending  the  lower  limit  of  the  above  integral  to 
negative  infinity  (the  lower  limit  already  spans  a  3^x  domain  of  the 
leading  edge  of  the  cloud).  The  integral  is  then  one  over  symmetric 
limits,  and  the  terms  involving  g^(t)  vanish.  Thus,  with  this 
approximation,  the  expression  for  the  unit  time  reference  dose  rate 
becomes 

Di(x,y)  =  kYFff(y,ta)g(ta)/V  (37) 

which  the  preceeding  derivation  shows  is  accurate  to  within  a  linear 
(not  constant)  approximation  for  g(t). 
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IV.  Validation  of  the  AFIT  Model  Activity  Deposition  Rate 

Obtaining  g(  t)  from  DELFIC  Data 

A  major  advantage  of  the  smearing  model  described  above  is  that 
calculation  of  fallout  pattern  dose  rate  contours  is  based  on 
independent  functions  describing  cloud  geometry  and  wind.  Thus,  in 
verifying  the  mathematics  of  the  smearing  model,  one  need  not  be 
concerned  with  its  ability  to  reproduce  either  observed  fallout 
patterns  or  patterns  calculated  by  some  comparison  standard  like 
DELFIC.  Since  the  current  formulation  of  the  AFIT  model  does  not  allow 
for  variable  winds,  good  agreement  with  actual  or  DELFIC  calculated 
patterns  could  be  expected  only  for  actual  winds  that  matched  the 
effective  wind  conditions  permitted  by  the  model.  Thus,  in  validating 
the  AFIT  model,  only  the  calculated  rate  of  arrival  of  activity  on  the 
ground,  g(t),  will  be  compared  with  the  prediction  of  a  reliable 
standard.  DELFIC' s  prediction  of  g(t)  is  chosen  as  that  standard. 

Given  a  well  formulated  g(t),  the  model's  ability  to  make 
reasonable  predictions  of  far  field  fallout  will  rest  in  its  treatment 
of  winds.  While  the  current  formulation  is  limited  to  a  wind  treatment 
neither  more  nor  less  restrictive  than  that  of  WSEG,  unlike  WSEG,  the 
AFIT  model  does  have  a  potential  for  incorporation  of  variable  wind 
effects.  Such  an  extension  of  the  model,  which  is  capable  of  producing 
a  curved  "hotline",  has  been  reported  by  Hopkins  (20). 

Although  DELFIC  was  not  originally  written  to  provide  a 
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calculation  of  g(t),  the  required  data,  in  the  form  of  cloud  wafer 
grounding  statistics  and  particle  activity  assignments,  are  available 
from  the  output  of  the  transport  and  particle  activity  modules.  The 
general  procedure  is  to  calculate  the  integral  of  g(t)  from  zero  time 
(t=0)  to  arrival  time  (t=t  ),  by  summing  the  contributions  of 
individual  cloud  wafers.  The  integral  curve  is  then  differentiated 
numerically  to  obtain  g(t). 


Comparison  of  AFIT  Model  and  DELFIC  g( t)  Results 

The  AFIT  smearing  model  is  validated  here  by  comparison  of  its 
predictions  for  the  ground  arrival  rate  of  fallout,  g(t),  with  those  of 
the  AFIT  version  of  DELFIC  (the  g(t)  calculation  was  never  incorporated 
into  the  original  ASA  version  of  the  code).  The  DELFIC  calculations 
assumed  Pu  fuel  fissioned  by  high-energy  (14  MeV)  neutrons;  the 
atmosphere  was  the  U.S.  Standard  Atmosphere  for  mid-latitudes,  spring 
or  fall.  Four  yields:  1,  10,  100,  and  1000  KT,  and  six  different 
lognormal  particle  size  frequency  distributions  were  selected;  100 
particle  size  classes  were  defined  by  DELFIC  for  each  case.  In 
addition  to  the  DELFIC  default  particle  size  distribution,  the  selected 
distributions  include  the  five  distributions  introduced  in  Chapter  II: 
the  four  previously  used  in  a  DELFIC  study  by  Norment  (21:25),  and  one 
developed  for  this  study  to  reproduce  the  RAND  (WSEG-10)  size-activity 
distribution.  Figure  IV-1  displays  the  two-component  size-activity 
distributions  calculated  for  the  AFIT  model  and  summarizes  the 
size-frequency  distribution  parameters  used  both  in  the  AFIT 
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Figure  IV-1  Size-Activity  Distributions  for  the  AFIT  Model 


calculation  and  by  DELFIC.  Figures  III-3  through  III-8  previously 
compared  the  AFIT  and  DELFIC  distributions. 

Figures  IV-2  (presented  here)  and  C-1  through  C-23  (Appendix  C) 
present  comparisons  of  the  DELFIC  results  for  the  activity  arrival  rate 
with  results  of  the  "pancake"  cloud  version  of  the  AFIT  smearing  model. 
Figures  IV-2  and  C-1  through  C-5  are  1000  KT  cases;  Figures  C-6  through 
C— 1 1  are  100  KT  cases;  Figures  C-12  through  C-17  are  10  KT  cases;  and 
Figures  C-18  through  C-23  are  1  KT  cases.  Several  features  of  these 
comparison  figures  are  immediately  evident.  First,  and  perhaps  most 
significant,  unlike  the  WSEG  model  calculation  presented  in  Figure 
III-1 ,  the  AFIT  smearing  model  closely  parallels  the  behavior  predicted 
by  DELFIC  for  g(t)  for  every  size-activity  distribution  at  every  yield; 
the  agreement  is  particularly  good  for  the  DELFIC  Default  and  Low  Yield 
Standard  cases.  The  trend  is  for  better  agreement  for  higher  yields 
and  later  times.  This  is  what  one  might  expect  in  view  of  the  AFIT 
model's  initial  condition  of  a  stabilized  cloud  over  ground  zero.  At 
higher  yields,  the  stem,  which  the  AFIT  model  lacks,  is  a  less 
significant  source  of  activity;  and  at  later  times,  one  is  concerned, 
in  both  models,  with  the  relatively  small  particles  which  began  their 
fall  from  the  cloud  cap.  The  trend  of  the  AFIT  model  to  predict 
slightly  higher  late  time  deposition  rates  than  DELFIC  is  explained  by 
the  AFIT  model's  assignment  of  a  somewhat  larger  fraction  of  activity 
to  the  smallest  particles. 

The  relative  significance  of  the  stem  depends  on  the  size-activity 
distribution,  as  well  as  on  the  yield.  The  worst  agreement  between 
DELFIC  and  the  AFIT  smearing  model  occurs  for  the  extreme  NRDL 
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Figure  IV-2  Deposition  Rates  for  the  DEFAULT  Distribution 
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distributions,  N61  and  C6l.  Since  these  two  distributions  have  rather 
large  geometric  deviations  (large  logarithmic  slopes),  they  include  a 
relative  abundance  of  very  large  particles.  In  the  DELFIC  simulation, 
such  particles  reside  primarily  in  the  stem,  from  which  they  are 
deposited  earlier  than  they  can  be  from  the  AFIT  model's  stabilized 
cloud.  The  net  effect  is  that  the  peak  in  the  AFIT  g(t)  is  generally 
higher  and  occurs  later  than  does  the  peak  in  the  DELFIC  g(t)  for  these 
distributions. 

At  early  times  (less  than  about  5  hours)  and  especially  at  low 
yields,  the  DELFIC  curves  display  slope  discontinuities  which  are  not 
present  in  the  curves  for  the  AFIT  model.  This  structure  arises  from 
the  numerical  differentiation  of  DELFIC' s  histogram  for  the  fraction  of 
total  activity  down  as  a  function  of  time.  The  histogram  cannot  be 
smooth  because  of  the  cloud  wafer  quadrature.  This  quadrature  also 
places  a  limit  on  the  ability  of  the  differentiation  to  resolve  the 
peak  in  the  arrival  rate  function.  DELFIC' s  quadrature  is  similarly 
responsible  for  the  occasional  late-time  discontinuities  in  the  arrival 
rate  slope. 

The  truncation  of  the  DELFIC  curves  at  about  70  hours,  or  earlier 
in  the  case  of  the  WSEG  or  HYSTND  distributions,  is  an  artifact  caused 
by  the  numerical  differentiation.  Since  the  arrival  rate  at  each  time 
was  taken  to  be  the  slope  of  the  deposited  activity  fraction  between 
that  time  and  the  next  interval,  and  only  points  which  fell  within  the 
zero  to  100  hour  plotting  interval  were  considered,  the  slope  was  never 
calculated  beyond  70  hours.  The  WSEG  and  the  HYSTND  cases  terminate 
earlier  because  the  deposition  rate  for  the  succeeding  time  point  was 
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zero. 


The  impact  of  the  differences  in  activity  grounding  rates  on  the 
total  fraction  of  activity  deposited  is  illustrated  by  Figures  IV-3, 
IV-4,  and  IV-5,  for  a  1  KT  yield,  with  the  DEFAULT,  WSEG-10,  and 
NRDLN61  distributions.  These  figures  demonstrate  that  the  most 
significant  differences  occur  prior  to  about  2  hours;  after  that,  the 
AFIT  and  DELFIC  deposition  fraction  curves  are  nearly  parallel  or 
slightly  convergent.  For  higher  yields  the  agreement  is  better; 
results  for  the  other  distributions  are  similar. 

To  the  extent  indicated  by  the  comparison  with  DELFIC  results,  the 
AFIT  model's  formulation  for  the  ground  arrival  rate  of  activity,  g(t), 
now  stands  validated.  Although  its  impact  on  g(t)  appears  to  be  minor, 
the  data  presented  above  suggest  that  a  significant  area  for 
improvement  in  the  g(t)  determination  is  in  the  modeling  of  cloud 
geometry. 
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Figure  IV-3  Activity  Deposition  for  the  DEFAULT  Distribution 


Figure  IV-H  Activity  Deposition  for  the  WSEG-10  Distribution 


Figure  IV-5  Activity  Deposition  for  the  NRDLN61  Distribution 


Sensitivity  of  the  AFIT  Model  to  Parameters 


of  the  Particle  Size-Activity  Distribution 

A  key  feature  of  the  cloud  smearing  approach  to  fallout  modeling 
is  the  use  of  independent  functions  to  describe  the  rate  of  arrival  of 
cloud  activity  on  the  ground  and  the  horizontal  spatial  distribution  of 
that  activity  in  the  cloud.  Together,  these  functions  produce  the 
fallout  footprint.  By  itself,  the  horizontal  cloud  distribution 
function  reflects  the  effects  of  wind,  wind  shear,  and  turbulence  on 
the  advective  transport  of  the  debris  cloud  as  it  drifts  downwind  from 
the  detonation  site.  On  the  other  hand,  the  activity  arrival  rate 
function,  g(t),  expresses  the  vertical  processes  involved  in  fallout. 
It  reflects  the  combined  relationship  between  the  initial  vertical 
distribution  of  activity,  h(z),  the  particle  size  distribution  of 
activity,  A(r),  and  the  physical  laws  describing  the  gravity  settling 
rates  of  fallout  particles.  In  the  smearing  model,  gravity  settling  is 
expressed  as  dr/dt:  the  initial  altitude  and  time-dependent  rate  of 

change  in  the  size  of  arriving  particles.  Thus  the  smearing  model 
separates  the  horizontal  process  of  advective  transport  of  fallout  from 
the  vertical  process  of  particle  settling.  This  feature  permits  the 
activity  arrival  function  to  be  used  as  a  measure  of  the  sensitivity  of 
the  model,  without  the  the  need  to  consider  the  complicating  effects  of 
horizontal  transport  under  the  influence  of  prevailing  winds. 

Among  the  three  factors  determining  the  activity  arrival  function, 


attention  is  focused  here  on  sensitivity  to  the  parameters  defining  the 


particle  size-activity  distribution.  The  following  observations  are 
offered  in  partial  justification  of  this  focus.  Since  the  particle 
fall  dynamics  incorporated  in  the  AFIT  model’s  evaluation  of  dr/dt,  or 
minor  modifications  thereof  (1:25),  are  considered  to  represent  the 
state-of-the-art  in  modeling  the  free-fall  of  spheres  in  a  viscous 
medium,  there  is  little  incentive  for  considering  older  approximations. 
Norment  (2 1:32)  has  previously  examined  the  impact  of  the  assumption 
that  fallout  particles  are  spherical,  rather  than  some  less  regular 
shape,  and  found  at  most  minor  impact  on  fallout  prediction.  With 
regard  to  the  vertical  distribution  of  activity,  it  was  observed 
previously  in  Figure  III-9  that  the  details  of  a  symmetric  distribution 
of  activity  about  the  cloud  center  height  have  a  negligible  effect  on 
the  ground  arrival  rate  of  activity,  g(t).  The  average  cloud  height, 
on  the  other  hand,  is  a  direct  measure  of  the  yield.  For  a  given 
size-activity  distribution,  it  determines,  through  the  dr/dt  relation, 
the  nature  of  the  arrival  rate  curve.  The  figure  just  cited  shows  that 
the  fractional  arrival  rate  of  activity  is  moderately  sensitive  to 
yield.  While  the  yield  varies  over  three  orders  of  magnitude,  the 
fractional  arrival  rate  varies  by  less  than  one  order  of  magnitude. 
Increasing  the  yield  decreases  the  magnitude  of  the  peak  in  the  gif' 
curve,  shifting  the  peak  to  slightly  later  times.  In  compensation,  the 
late  time  arrival  rates  are  higher  for  higher  yields. 

Although  the  distribution  of  activity  with  fallout  particle  sice 
is  of  prime  significance  in  the  calculation  of  g(t),  determination  of 
an  appropriate  distribution  is  uncertain.  Experimental  support  for 
distributions  such  as  the  DELFIC  default  particle  size-frequency 


distribution  is  generally  meager  and  potentially  biased  by  sampling 
techniques.  The  DELFIC  default  distribution  is  based  on  Nevada  test 
site  debris.  Its  applicability  to  predicting  fallout  resulting  from  a 
detonation  on  hard  rock  or  in  a  large  city  is  dubious.  In  short,  one 
of  the  most  important  factors  in  the  prediction  of  fallout  is  largely 
unknowable.  Thus,  it  is  important  to  understand  the  impact  of  the 
uncertainties  in  distribution  of  activity  with  particle  size  on  fallout 
prediction.  This  is  impossible  with  the  current  mainstay  of 
operational  type  studies,  WSEG-10,  because  in  WSEG  the  activity 
distribution  and  fall  rates  have  been  pre-combined.  With  the  AFIT 
smearing  model,  however,  the  impact  of  varying  the  individual 
parameters  is  directly  observable. 

The  sensitivity  of  the  AFIT  smearing  model  to  the  distribution  of 
activity  with  particle  size  was  examined  by  observing  the  effects  on 
activity  grounding  times  of  varying  the  distribution  parameters.  The 
extreme  sensitivity  of  the  model  to  such  variations  is  readily  apparent 
from  an  examination  of  the  deposition  rate  data  presented  in  Figures 
IV-2  and  C-1  through  C-23  as  a  part  of  the  model  validation.  Here, 
with  the  DELFIC  default  size-frequency  distribution  as  a  baseline,  each 
of  the  three  parameters  that  define  the  AFIT  model's  two-component 
size-activity  distribution  are  varied:  the  refractory  activity 
fraction,  that  is,  the  fraction  of  activity  distributed  volumetrically , 
f  ,  the  particle  size-frequency  median  diameter,  do,  and  the  factor 
defining  the  dispersion  of  particle  sizes  about  that  median,  the 
geometric  standard  deviation,  s. 

Figure  III-2  showed  that  values  of  the  .’efractory  activity 
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fraction,  f^,  calculated  on  the  basis  of  DELFIC's  particle  activity 
assignment  data,  were  insensitive  to  yield  and  only  slightly  dependent 
on  choice  of  particle  size  distribution.  Although  the  observed  range 
of  f  was  only  from  0.66  to  0.7^,  the  DELFIC  calculations  all  assumed 
fission  of  Pu  fuel  by  14  MeV  neutrons.  Different  fuels,  different 
neutron  energy  spectra,  or  the  presence  of  significant  induced  activity 
species  could  significantly  alter  the  fractionation  behavior  of  the 
debris.  Here  the  entire  range  of  possible  values  for  the  refractory 
activity  fraction  is  considered.  Figure  V-1  displays  the  effects  on 
the  size-activity  distribution  of  varying  the  fractionation  parameter, 
f  ,  from  0.0,  corresponding  to  a  totally  surface  distribution  of 
activity,  to  1.0,  corresponding  to  a  totally  volumetric  (refractory) 
distribution.  Although  this  leads  to  a  variation  in  the  size-activity 
medians  spanning  the  range  from  surface  to  mass  distribution  median, 
the  dispersion  about  those  medians  is  similar  to  that  of  the 
size-frequency  distribution  itself.  Increasing  the  refractory  activity 
fraction  has  the  effect  of  increasing  the  median  of  the  size-activity 
distribution.  Figures  V-2  through  V-5  show  the  effects  of  these 
variations  on  the  fractional  arrival  rate  of  activity  for  yields 
ranging  from  1  to  1000  KT.  For  each  yield,  the  variations  in  the 
refractory  activity  fraction  produce  a  family  of  curves  lying  within  an 
envelope  bounded  by  the  surface  and  volume  distribution  results.  At 
early  times,  the  more  refractory  distributions  lead  to  higher  peak 
deposition  rates;  at  late  times  they  lead  to  lower  deposition  rates. 
The  spread  in  peak  rates  is  about  one  order  of  magnitude.  The 


cross-over  from  "early"  to  "late",  the  nodal  point  in  the  envelope, 


Figure  V-1  Activity  Distribution  Sensitivity  to  Fractionation 
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Figure  V-3  Deposition  Rate  Sensitivity  to  Fractionation 
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occurs  at  the  time  corresponding  to  deposition  of  particles  whose  size 
is  the  geometric  mean  of  the  surface  and  volume  distribution  medians, 
the  median  of  the  2.5th  moment  of  the  size-frequency  distribution.  To 
either  side  of  the  nodal  time,  the  variation  of  g(t)  with  f  is 
monotonic.  The  node  itself  shifts  to  later  times  as  yield  increases 
and  the  time  required  for  any  fixed  size  particle  to  be  deposited 
increases.  Thus,  the  magnitude  and  timing  of  the  peak,  as  well  as  the 
late  time  slope  of  the  activity  arrival  rate  function  are  all 
influenced  by  the  choice  of  the  refractory  activity  fraction. 

Figure  V-6  displays  the  effects  on  the  size-activity  distribution 
of  variations  in  the  size-frequency  median  diameter,  do,  from  one  tenth 
to  ten  times  the  DELFIC  default  value  of  0.407  ym.  This  leads  to  a 
comparable  two  order  of  magnitude  span  between  the  medians  of  the 
corresponding  size-activity  distributions.  Again,  the  dispersion  of 
activity  about  the  size-activity  medians  is  similar  to  that  of  particle 
sizes  about  the  size-frequency  median.  Since  the  impact  of  these 
variations  is  a  more  extreme  example  of  the  same  effects  as  those 
produced  by  fractionation  variations,  but  without  the  surface  and 
volume  distribution  bounds,  only  the  deposition  rate  results  for  the 
1000  KT  example  are  presented  (Figure  V-7).  As  with  the  previous 
fractionation  variations,  there  is  approximately  a  one-to-one 
relationship  between  peak  arrival  rate  differences  and  size-activity 
median  variations.  The  peak  in  the  activity  arrival  rate  curve  for  the 
distribution  with  the  largest  median  occurs  earlier  and  is  higher  and 
steeper  than  the  peaks  for  the  other  distributions.  Since  this 
distribution  has  relatively  little  activity  associated  with  the  smaller 
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Figure  V-6  Activity  Distribution  Sensitivity  to  Size  Median 


particles,  its  activity  arrival  rate  declines  more  rapidly  at  times 
beyond  the  peak  than  does  the  rate  for  the  other  distributions.  On  the 
other  hand,  the  arrival  rate  curve  for  the  distribution  with  the 
smallest  median  rises  to  a  broad  maximum.  Since  this  distribution 
contains  an  abundance  of  activity  on  small  particles  with  long  fall 
times,  the  activity  arrival  rate  decreases  slowly  with  time. 

Figure  V-8  displays  the  effects  on  the  size-activity  distribution 
of  variations  in  the  dispersion  of  particle  sizes  about  the  median,  as 
measured  by  the  geometric  standard  deviation  of  the  distribution,  s. 
Here,  s  has  been  allowed  to  vary  by  one  unit  above  and  below  the 
standard  value  of  4.0  for  the  DELFIC  default  distribution.  In 
appreciating  the  effect  (nearly  an  order  of  magnitude)  that  this  small 
variation  in  the  geometric  standard  deviation  has  on  the  size-activity 
median,  it  is  helpful  to  recall  the  expression  for  the  medians  of  the 
various  moment  distributions  of  the  size-frequency  distribution, 
equation  (32).  Although  the  AFIT  model's  two-component  size-activity 
distribution  is  not  lognormal,  but  rather  the  sum  of  two  lognormal 
distributions,  with  a  refractory  fraction  of  0.68,  its  median  will  lie 
between  that  of  the  2.5th  and  third  moment  distributions.  Since  the 
median  diameter  of  the  size-frequency  distribution  is  0.407,  even  the 
median  of  the  surface  distribution  is  dominated  by  the  second  term  in 
equation  (32)  involving  the  geometric  standard  deviation.  In  fact,  the 
size-activity  distribution  with  the  small  dispersion  of  particle  sizes 
(s=3)  has  a  median  smaller  than  that  of  the  distribution  of  surface 
area  for  the  baseline  distribution,  and  the  distribution  with  the  large 
dispersion  of  particle  sizes  (s=5)  has  a  median  larger  than  that  of  the 
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Figure  V-8  Activity  Distribution  Sensitivity  to  Size  Dispersion 


distribution  of  volume  for  the  baseline.  As  a  result,  Figure  V-9  shows 
that  the  small  and  large  deviation  deposition  rate  curves  bracket  the 
results  presented  previously  in  Figure  V-2,  which  considered 
fractionation  variations. 

This  analysis  has  shown  that  the  choice  of  the  particle 
size-frequency  distribution  and  the  partition  of  activity  between 
volume  and  surface  modes  of  distribution  are  of  paramount  significance 
to  fallout  model  predictions. 
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VI .  Conclusion 


This  investigation  has  pursued  two  methods  for  improving  far  field 
fallout  prediction  capability,  identified  above  as  cloud  wafer 
smoothing  and  whole  cloud  smearing.  In  the  former  approach,  the 
advection-plus-settling  transport  algorithm  in  the  state-of-the-art 
disc  tosser  code,  DELFIC,  as  received  from  Atmospheric  Sciences 
Associates,  was  upgraded  to  improve  overlap  of  grounded  cloud  wafers  by 
spreading  their  contents  over  a  ground  area  that  more  accurately 
represented  the  true  dispersal  of  the  wafer  contents.  By  maintaining 
the  contiguous  nature  of  cloud  wafer  boundaries  from  cloud  rise  through 
ground  deposition,  it  was  possible  to  guarantee  that  no  voids  would 
exist  in  downwind  fallout  patterns,  regardless  of  the  distance  from  the 
burst  location.  The  approach  is  only  partially  successful,  however,  in 
that  realism  in  far  field  predictions  by  a  disc  tosser  requires  the 
overlap  at  each  ground  point  of  several  cloud  wafers.  Otherwise,  the 
far  downwind  dose  rate  contour  shapes  are  determined  primarily  by  the 
assumed  horizontal  distribution  of  activity  in  each  wafer  of  the 
stabilized  cloud. 

In  the  second  method  of  fallout  modeling,  cloud  wafer  smoothing  is 
carried  to  the  limit  by  treating  the  whole  cloud  as  a  single  wafer 
composed  of  activity  distributed  in  a  continuous  manner  over  the 
'~arrier  soil  particle  distribution.  A  technique  was  developed  to  smear 
'-loud  contents  continuously  over  the  ground.  Although  the  smearing 


model  retains  the  potential  for  incorporation  of  a  vertical  wind 
profile,  within  the  scope  of  this  investigation,  that  capability  was 
sacrificed,  with  attention  focused  on  calculation  of  the  fractional 
arrival  rate  of  activity  on  the  ground,  g(t),  by  a  method  not 
restricted  as  disc  tossers  are,  to  regions  relatively  close  to  ground 
zero.  The  method  of  calculating  g(t)  described  here  and  dubbed  the 
AFIT  smearing  model  has  been  shown  to  produce  g(t)  values  which  agree 
well  with  DELFIC  predictions. 

By  providing  an  accounting  for  fractionation  and  a  means  to 
perform  sensitivity  analyses,  the  AFIT  model  introduces  capabilities 
not  included  in  earlier  smearing  models.  The  model  achieves  a  distinct 
advantage  by  decoupling  the  effects  of  the  various  input  parameters. 
The  fractional  ground  arrival  rate  of  activity,  g(t),  is  determined  by 
the  initial  1  distribution  of  activity,  by  the  size-frequency 
distributior  ’ained  soil  particles,  by  fractionation  of  activity 
between  volume  .id  surface  distributions,  and  by  particle  settling 
rates.  In  contrast,  g(t)  is  completely  independent  of  the  horizontal 
distribution  of  activity  in  the  cloud,  of  horizontal  turbulent  growth 
of  the  cloud,  and  of  wind  speed  and  shear.  As  a  result,  since  the 
parameters  which  influence  g(t)  are  directly  and  individually 
accessible,  the  AFIT  model  becomes  a  useful  tool  in  surveying  the 
sensitivity  of  predictions  to  variations  in  those  parameters. 

The  sensitivities  of  the  AFIT  smearing  model  to  variations  in 
yield  (vertical  distributin'  of  activity)  and  distribution  of  activity 
with  particle  size  (including  fractionation  effects)  were  implicitly 


demonstrated  by  the  various  comparisons  with  DELFIC  results  during  the 
validation  phase;  sensitivity  was  explicitly  considered  in  Chapter  V, 
through  variation  of  the  parameters  of  the  two-component  size-activity 
distribution:  the  fractionation  parameter,  the  size-frequency  median, 
and  the  size-frequency  geometric  standard  deviation.  The  peak  time  and 
maximum  value,  as  well  as  the  late  time  slope  of  the  fractional 
activity  arrival  rate  curve  were  found  to  be  extremely  sensitive  to  the 
geometric  standard  deviation  of  the  size-frequency  distribution.  This 
parameter  directly  affects  both  the  median  and  dispersion  of  the  size- 
activity  distribution.  When  the  size-frequency  median  alone  or  the 
refractory  activity  fraction  are  varied,  the  median  of  the  resulting 
activity  distribution  is  altered  proportionately,  while  the  dispersion 
of  the  activity  distribution  remains  similar  to  the  dispersion  of  the 
frequency  distribution.  As  the  refractory  activity  fraction  varies 
between  zero  and  one,  the  size-activity  median  and  the  peak  magnitude 
of  the  fractional  arrival  rate  vary  by  an  order  of  magnitude  or  more. 
Thus,  the  AFIT  model  is  extremely  sensitive  to  the  particle 
size-frequency  spectrum  and  to  the  partition  of  activity  between  volume 
and  surface  modes  of  distribution. 

The  validation  of  the  AFIT  model  prediction  for  g(t)  by  comparison 
to  DELFIC  results  demonstrates  the  model’s  competence  as  a  method  for 
predicting  the  ground  arrival  rate  of  activity.  Since  it  is  capable  of 
treating  fallout  transport  without  resort  to  numerical  quadrature,  it 
can  suffer  none  of  the  artificial  structure  problems  inherent  in  disc 
tosser  codes.  In  its  treatment  of  winds,  however,  it  is  neither  more 


vv 


nor  less  competent  than  its  predecessor,  WSEG.  This  is  a  productive 

area  for  further  model  improvement.  In  view  of  the  strong  sensitxvity 

to  the  fractionation  parameter,  a  direct,  first-principles  calculation 

of  f  would  also  be  a  welcome  addition  to  the  model, 
r 
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Appendix  A 

Fallout  Patterns  Produced  by  the  AFIT  Version  of  DELFIC 

The  following  figures  (A-1  through  A-24)  present  fallout  contour 

pattern  predictions  generated  by  the  AFIT  version  of  DELFIC.  The 

particle  size-frequency  distribution  parameters  were  presented  in  Table 

II-1.  Figure  II-4  showed  the  corresponding  DELFIC  activity 

239 

assignments.  The  activity  calculations  assumed  fission  of  Pu  by  14 
MeV  neutrons  (Fission  ID  =  P239HE),  and  a  3oil  particle  melting  point 
of  2200°K.  The  quantities  wD(MED)n  and  "S(GEO),r  are  parameters  of  a 
lognormal  least  squares  fit  to  the  resulting  DELFIC  size-activity 
assignments.  The  atmospheric  properties  used  were  those  of  the  U.  S. 
Standard  Atmosphere  for  spring  or  fall  at  mid-latitudes  (Atmosphere  ID 
=  MIDLSF ) .  The  winds  used  were  those  measured  at  Washington,  DC,  on 
October  3,  1968  (Table  II-2).  Each  figure  is  identified  by  a  code 
derived  from  the  wind  profile,  yield,  and  size  distribution.  The  wind 
information  shows  date  and  location;  "03W"  means  October  3,  1968,  at 
Washington,  DC.  Yields  are  indicated  in  kilotons,  with  a  trailing  "S" 
for  surface  burst.  The  DELFIC  default  size  distribution  is  "DEF" ;  the 
WSEG-10,  SIDAC,  or  RAND  is  "WSEG";  the  others  are  identified  on  the 
basis  of  the  particle  size-frequency  median,  for  example,  "22P7"  for 
22.7  Wm.  Thus,  "03W1000S22P7",  indicates  a  1000  KT  surface  burst 
calculation,  using  winds  observed  on  October  3,  1968,  at  Washington, 
DC,  and  employing  the  size  distribution  having  a  median  of  22.7  Wm. 


Figure  A-1  Fallout  Contours  for  O3WI000SDEF 


Figure  A-2  Fallout  Contours  for  O3WI00QSP3I4 


Figure  A-3  Fallout  Contours  for  03Wt 


Figure  A-4  Fallout  Contour*  for  03W1000SWSBC 


Figure  A -5  Fallout  Contours  for  03WIOOOSPOB49 


Figure  A-6  Fallout  Contours  for  O3WI000SI 
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figure  A-7  Fallout  Contours  for  O3WI00SDEX 
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Figure  A-9  Fallout  Contours  for  03W100S22P? 
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Figure  A-10  Fallout  Contours  for  O3WI00SWSBC 
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Figure  A-13  Fallout  Contours  for  O3WI0SDEF 


figure  A- 14  fallout  Contours  for  03WI0SP3I4 
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Figure  A-l6  Fallout  Contour*  for  03WIOSWSBG 


Figure  A- 17  Fallout  Contour*  for  O3WI0SP0849 
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DISTANCE  EAST  of  C-Z  (meters) 
-18  Fallout  Contours  for  O3WI0SPO 
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DISTANCE  EAST  of  C-Z  (meters) 
Figure  A-19  Fallout  Contours  for  03WISDEF 


Figure  A-20  Fa f foul  Contours  for  03WI5P3I4 


DISTANCE  EAST  of  G-Z  (meters) 

A -21  Fallout  Contours  for  03WISZ2P7 


Figure  A-22  Fallout  Contours  for  03WISWSBC 
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Figure  A-23  Fallout  Contours  for  03WISP0849 


DISTANCE  EAST  of  C-Z  (meters) 

Figure  A-24  Fallout  Contours  for  O3WISP0C 
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Deposition  Rates  for  the  HYSTND  Distribution 
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Figure  C-3  Deposition  Rates  for  the  WSEG-10  Distribution 
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Figure  C-4  Deposition  Rates  for  the  NRDLC61  Distribution 
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Figure  C-5  Deposition  Rates  for  the  NRDLN61  Distribution 
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Figure  C-7  Deposition  Rates  for  the  LYSTND  Distribution 
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Figure  C-8  Deposition  Rates  for  the  HYSTND  Distribution 
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Figure  C-9  Deposition  Rates  for  the  WSEG-10  Distribution 
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Figure  C-10  Deposition  Rates  for  the  NRDLC61  Distribution 
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Figure  C-13  Deposition  Rates  for  the  LYSTND  Distribution 
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Figure  C-17  Deposition  Rates  for  the  NRDLN61  Distribution 


■=ramc  T-  1  kt,  vo-'tOT 

io«*m  Y-  1  KT.  d.-o.4Q7 


'I . I . I . . . . . . 

aa  40  so  ao  to  ao  oo 


TIME  (hours) 


Figure  C-18  Deposition  Rates  for  the  DEFAULT  Distribution 


Figure  C-19  Deposition  Rates  for  the  LYSTND  Distribution 
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